#### FIGURE E.2: NEGATIVE CONSEQUENCES PERSISTENCE
#### Study 3 waves for welfare, crime, and employment outcomes

rm(list = ls())
source("./2_code/00_setup.R")

#### STUDY 3, WAVE 1 ####

data3 <- fread(paste0(data_path, "data_study3.csv"), header = TRUE)
data3$employed <- as.factor(data3$employed)

# Welfare
lin_welfare_economic3 <- lm_lin(immigrant_welfare ~ treat_economic, covariates = ~
                                  edad + employed + male, data = data3[(treat=='economic' | treat=='control')])
lin_welfare_humanitarian3 <- lm_lin(immigrant_welfare ~ treat_humanitarian, covariates = ~
                                      edad + employed + male, data = data3[(treat=='humanitarian' | treat=='control')])

# Crime
lin_crime_economic3 <- lm_lin(immigrant_crime ~ treat_economic, covariates = ~
                                edad + employed + male, data = data3[(treat=='economic' | treat=='control')])
lin_crime_humanitarian3 <- lm_lin(immigrant_crime ~ treat_humanitarian, covariates = ~
                                    edad + employed + male, data = data3[(treat=='humanitarian' | treat=='control')])

# Employment
lin_employ_economic3 <- lm_lin(immigrant_employment ~ treat_economic, covariates = ~
                                 edad + employed + male, data = data3[(treat=='economic' | treat=='control')])
lin_employ_humanitarian3 <- lm_lin(immigrant_employment ~ treat_humanitarian, covariates = ~
                                     edad + employed + male, data = data3[(treat=='humanitarian' | treat=='control')])


#### STUDY 3, WAVE 2 ####

data3b <- fread(paste0(data_path, "data_study3_wave2.csv"), header = TRUE)

# Welfare
lin_welfare_economic3b <- lm_lin(immigrant_welfare ~ treat_economic, covariates = ~
                                   edad + employed + male, data = data3b[(treat=='economic' | treat=='control')])
lin_welfare_humanitarian3b <- lm_lin(immigrant_welfare ~ treat_humanitarian, covariates = ~
                                       edad + employed + male, data = data3b[(treat=='humanitarian' | treat=='control')])

# Crime
lin_crime_economic3b <- lm_lin(immigrant_crime ~ treat_economic, covariates = ~
                                 edad + employed + male, data = data3b[(treat=='economic' | treat=='control')])
lin_crime_humanitarian3b <- lm_lin(immigrant_crime ~ treat_humanitarian, covariates = ~
                                     edad + employed + male, data = data3b[(treat=='humanitarian' | treat=='control')])

# Employment
lin_employ_economic3b <- lm_lin(immigrant_employment ~ treat_economic, covariates = ~
                                  edad + employed + male, data = data3b[(treat=='economic' | treat=='control')])
lin_employ_humanitarian3b <- lm_lin(immigrant_employment ~ treat_humanitarian, covariates = ~
                                      edad + employed + male, data = data3b[(treat=='humanitarian' | treat=='control')])


#### PREPARE DATA FOR PLOTTING ####

types_economic <- data.frame(
  treatment = rep('economic', 2),
  Study = c('Wave 1', 'Wave 2')
)

types_humanitarian <- data.frame(
  treatment = rep('humanitarian', 2),
  Study = c('Wave 1', 'Wave 2')
)


### Economic treatment, welfare ###

out_economic_welfare <- rbind(
  tidy(lin_welfare_economic3) %>% filter(term=='treat_economic') %>% dplyr::select(estimate, std.error),
  tidy(lin_welfare_economic3b) %>% filter(term=='treat_economic') %>% dplyr::select(estimate, std.error)
)

out_economic_welfare <- cbind(out_economic_welfare, types_economic)


### Humanitarian treatment, welfare ###

out_humanitarian_welfare <- rbind(
  tidy(lin_welfare_humanitarian3) %>% filter(term=='treat_humanitarian') %>% dplyr::select(estimate, std.error),
  tidy(lin_welfare_humanitarian3b) %>% filter(term=='treat_humanitarian') %>% dplyr::select(estimate, std.error)
)

out_humanitarian_welfare <- cbind(out_humanitarian_welfare, types_humanitarian)


### Economic treatment, crime ###

out_economic_crime <- rbind(
  tidy(lin_crime_economic3) %>% filter(term=='treat_economic') %>% dplyr::select(estimate, std.error),
  tidy(lin_crime_economic3b) %>% filter(term=='treat_economic') %>% dplyr::select(estimate, std.error)
)

out_economic_crime <- cbind(out_economic_crime, types_economic)


### Humanitarian treatment, crime ###

out_humanitarian_crime <- rbind(
  tidy(lin_crime_humanitarian3) %>% filter(term=='treat_humanitarian') %>% dplyr::select(estimate, std.error),
  tidy(lin_crime_humanitarian3b) %>% filter(term=='treat_humanitarian') %>% dplyr::select(estimate, std.error)
)

out_humanitarian_crime <- cbind(out_humanitarian_crime, types_humanitarian)


### Economic treatment, employment ###

out_economic_employ <- rbind(
  tidy(lin_employ_economic3) %>% filter(term=='treat_economic') %>% dplyr::select(estimate, std.error),
  tidy(lin_employ_economic3b) %>% filter(term=='treat_economic') %>% dplyr::select(estimate, std.error)
)

out_economic_employ <- cbind(out_economic_employ, types_economic)


### Humanitarian treatment, employment ###

out_humanitarian_employ <- rbind(
  tidy(lin_employ_humanitarian3) %>% filter(term=='treat_humanitarian') %>% dplyr::select(estimate, std.error),
  tidy(lin_employ_humanitarian3b) %>% filter(term=='treat_humanitarian') %>% dplyr::select(estimate, std.error)
)

out_humanitarian_employ <- cbind(out_humanitarian_employ, types_humanitarian)


#### PLOTTING ####

fla <- rbind(
  out_economic_welfare,
  out_humanitarian_welfare,
  out_economic_crime,
  out_humanitarian_crime,
  out_economic_employ,
  out_humanitarian_employ
)

data_plot <- data.frame(
  fla,
  Outcome = c(rep("ATE Welfare Burden", 2),
              rep("ATE Welfare Burden", 2),
              rep("ATE Crime", 2),
              rep("ATE Crime", 2),
              rep("ATE Employment", 2),
              rep("ATE Employment", 2)),
  Treatment = c(rep("Economic Narrative", 2), rep("Humanitarian Narrative", 2),
                rep("Economic Narrative", 2), rep("Humanitarian Narrative", 2),
                rep("Economic Narrative", 2), rep("Humanitarian Narrative", 2))
)

# Reorder factors
data_plot$Study <- factor(data_plot$Study, levels = c("Wave 2", "Wave 1"))
data_plot$Outcome <- factor(data_plot$Outcome, levels = c("ATE Welfare Burden", "ATE Crime", "ATE Employment"))

plot_neg_waves <- ggplot(data_plot, aes(y = Outcome, x = estimate, color = Study, group = Study)) +
  geom_errorbarh(aes(xmin = estimate - 1.96 * std.error, xmax = estimate + 1.96 * std.error),
                 height = 0, position = position_dodge(width = 0.3), linewidth = 0.9) +
  geom_errorbarh(aes(xmin = estimate - 1.645 * std.error, xmax = estimate + 1.645 * std.error),
                 height = 0, position = position_dodge(width = 0.3), linewidth = 1.9) +
  geom_point(size = 3, position = position_dodge(width = 0.3)) +
  facet_grid(~ Treatment) +
  theme_minimal() +
  ylab("") + xlab("") +
  scale_x_continuous(limits = c(-0.3, 0.1)) +
  geom_vline(xintercept = 0, linetype = "dashed", color = "red") +
  guides(color = guide_legend(reverse = TRUE)) +
  labs(color = "Study 3") +
  theme(
    strip.text = element_text(size = 16),
    axis.text.y = element_text(size = 13),
    axis.text.x = element_text(size = 13),
    legend.text = element_text(size = 13),
    legend.title = element_text(size = 13)
  )

plot_neg_waves

ggsave(filename = paste0(plot_path, "figure_E2.png"), plot = plot_neg_waves, width = 10, height = 5, dpi = 600)